An advanced multipole model for (216) Kleopatra triple system

To interpret adaptive-optics observations of (216) Kleopatra, we need to describe an evolution multiple moons, orbiting extremely irregular body and including their mutual interactions. Such orbits are generally non-Keplerian orbital elements not constants. Consequently, use a modified $N$-body integrator, which was significantly extended include the multipole expansion gravitational field up order $\ell = 10$. Its convergence verified against `brute-force' algorithm. We computed coefficients $C_{\ell m},S_{\!\ell m}$ for Kleopatra's shape, assuming a~constant bulk density. For solar-system applications, it also necessary implement variable distance geometry observations. Our $\chi^2$ metric then accounts absolute astrometry, relative astrometry (2nd moon with respect 1st), angular velocities, silhouettes, constraining pole orientation. This allowed us derive two moons. Using both archival astrometric data new VLT/SPHERE (ESO LP 199.C-0074), were able identify true periods $P_1 (1.822359\pm0.004156)\,{\rm d}$, $P_2 (2.745820\pm0.004820)\,{\rm d}$. They orbit very close 3:2 mean-motion resonance, but osculating eccentricities too small compared other perturbations (multipole, mutual), so that regular librations critical argument present. The resulting mass $m_1 (1.49\pm0.16)\cdot10^{-12}\,M_\odot$ or $2.97\cdot10^{18}\,{\rm kg}$, is lower than previously thought. An implication explained in accompanying paper (Marchis et al.) Kleopatra critically rotating body.